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Abstract 

We propose an ensemble algorithm, which provides a new approach for evaluating and summing 
up a set of function samples. The proposed algorithm is not a quantum algorithm, insofar it does 
not involve quantum entanglement. The query complexity of the algorithm depends only on the 
scaling of the measurement sensitivity with the number of distinct spin sub-ensembles. From a 
practical point of view, the proposed algorithm may result in an exponential speedup, compared 
to known quantum and classical summing algorithms. However in general, this advantage exists 
only if the total number of function samples is below a threshold value which depends on the 
measurement sensitivity. 
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I. INTRODUCTION 



In this paper, we propose an ensemble algorithm for evaluating and summing an arbitrary 
function, as an alternative to the quantum algorithms that are currently believed to be the 
most efficient algorithms available, in terms of query complexity. 

The possibility of speeding up the evaluation and summing of a large number of function 
samples was first noted by Abrams and Williams , who suggested calculating numerical in- 
tegrals and stochastic processes using quantum algorithms. Quantum algorithms exploit the 
inherent parallelism offered by entangled quantum states, to perform certain computational 
tasks much more efficiently than classical devices, using either pure or pseudopure quantum 
states H, ^] that are tensor products of multiple qubits. The numerical value of an integral 
is evaluated by employing either the mean estimation algorithm devised by Grover || to 
calculate the mean of a discrete set of numbers, or by using the quantum counting algorithm 
proposed by Brassard, Hoyer, and Tapp || to determine the number of elements that fulfill 
a specified condition. Both of these approaches rely on a generalization of Grover's search 
algorithm, resulting in a quadratic speedup in comparison with classical randomized (Monte 
Carlo) algorithms, and an exponential speedup in comparison with classical deterministic 
algorithms for a single processor. 

A systematic comparison of optimal summation of finite sequences and continuous- 
function integration for deterministic, randomized, and quantum algorithms has been done 
by Heinrich and Novak j|, |5|. They have examined the query complexity of quantum integra- 
tion for different classes of integrand functions, assuming that the critical quantum speedup 
is obtained by using one of the two quantum summing algorithms mentioned above. 

Recently an alternative paradigm for computing has been suggested by Madi, Br- 
uschweiler, and Ernst, which operates on ensembles i.e., mixed states of identical spin sytems, 



using a spin Liouville space formalism [[12] , These types of ensemble algorithms are not 
quantum algorithms, insofar they do not involve entanglement of quantum states. Through- 
out this paper, "mixed states" describe a statistical ensemble, not individual systems, so 
that each element of the ensemble performs part of the computation, in the same way as a 
classical parallel computer. 

This new paradigm exploits the parallelism offered by simultaneously acting on linear 
combinations of many different input states in an ensemble of spins. Thus in general, 
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ensemble computing requires an exponentially larger set of memory resources to encode the 
same number of distinct input states compared to quantum computing with pure states. 

While ensemble computing requires more physical resources, it holds the important ad- 
vantage that it is insensitive to the decoherence time of the spins, which is an outstanding 
limiting factor for quantum computations involving entangled states. Moreover, ensemble 
algorithms can be exponentially faster than the equivalent quantum algorithms, for adequate 
measurement sensitivities, so a trade-off exists between memory and speed capabilities. 

We proceed to present a new approach to summing up function samples using an ensem- 
ble algorithm, and discuss its query complexity. We only consider the query complexity of 
the algorithm, i.e., the number of function evaluations that are performed, since the overall 
computational complexity will depend on the actual function that is being evaluated. At 
present, the most feasible physical implementation of this summing algorithm would rely on 
NMR technology, though any physical system of spins can be used in principle. The ideal 
physical system would allow us to have full control over a very large number of spins, in order 
to satisfy the large memory requirements. In the Discussion section, we comment on the 
application of the proposed summing algorithm to evaluating the mean of a continuous func- 
tion, and as a corollary, on estimating the definite integral of a continuous mult i- dimensional 
function. 

II. ENSEMBLE SUMMING ALGORITHM 
A. Statement of the Problem 

Let / : {1, 2, ... , N} — > [0, 1] be a real- valued function defined on a discrete set of samples 
comprised of N = 2 n points. The function / may be known analytically or it may be the 
result of an explicit or hidden numerical computation. The latter case is known as an oracle. 
We want to evaluate efficiently the sum Sn, 

N 

sw = E/(0- (i) 

i=i 

Here efficiency is understood in relation to the query complexity of the algorithm. Indeed, 
when N is large and the function evaluation is costly in terms of computational complexity, 
reducing the number of function evaluations is critical. 
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We assume that the algorithm is to be implemented in a physically realizable system 
consisting of a finite number of two-valued spins. To accomodate the N input values, we 
need n spins in the input register. The finiteness of the system and the discreteness of the 
spin states implies that we have to approximate the set of function values, {/(«)}, with a 
set of finite-precision values {/j G [0, 1]} for i — 1, 2, . . . , N. For the sake of simplicity we 
shall consistently use the same notation for f(i) and its fc-digit approximation. Since the 
meaning of the formulas is explained in the text, there should be no ambiguity. The number 
of spins k available in the output register will specify the minimal precision, 5 = 2~ k , for 
these values. Therefore we are actually evaluating the sum Sjv,fc, 

N 

SN,k = /»' ( 2 ) 

1=1 

which converges exponentially fast to the sum Sn, as we increase the number of spins k 
in the output register. Thus if we can evaluate <Sjv,fc efficiently, we can also evaluate Sn 
efficiently 

B. Outline of the Algorithm 

The proposed algorithm has three main steps. The first step consists of preparing an 
ensemble mixture of input states representing the numbers i = 1,2,..., AT. In the second 
step, the function / is applied to the input states, using a single transformation Uf to 
perform the function evaluation for every input state i at once. This parallel application 
results in an ensemble mixture which contains all of the values j\ in the output register. 
Finally, measurement of the output register automatically averages the contributions from 
the entire ensemble, yielding a signal proportional to the approximate sum, 5V,fc- 

Step 1 - Initialization 

We initialize the n-spin input register in an equally- weighted mixed state , 

1 N 

Pin = "jTfSl* >n< *\ n > ( 3 ) 
iV i=\ 

which accounts for all N = 2 n possible states. The mixed state is a density operator, 
which can be represented in spin Liouville space by a density matrix that has non-zero 
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elements only on its diagonal. The off-diagonal elements are all zero, indicating the absence 
of quantum coherence between any of the states \i > n . 

We note again that this mixed state describes a statistical ensemble, thus the initialization 
is equivalent to assigning each of the input values i = 1,2, ... ,N, to one of N classical 
processors in a parallel computer. 

For example, in an NMR implementation, the states \i > n correspond to the eigenstates 
of the Zeeman Hamiltonian created by a strong external magnetic field fl4| . The ket states, 
\i > n , can also be written in terms of individual spins, 

\i > n = | Oji > ®\a i2 > <g> ... <g> \a in > (4) 

where (an, a^, • • • , din) £ {0, 1} are the digits of the number (i — 1) in binary format. The 
bra states, < i\ n , are the dual of the ket states. The state |0 > denotes a spin "up" and the 
state |1 > denotes a spin "down". At room temperature, the thermal equilibrium state of 
an n-spin ensemble in an NMR experiment closely approximates the desired initial state p^ 
15]. The thermal state is equal to the sum of pj™ and a traceless deviation density matrix, 
with zero off-diagonal terms. The error introduced by using a thermal state instead of the 
equally-weighted mixed state is addressed in the Discussion. 

We also assume that we have available an output register with k spins, which is capable 
of encoding the real-numbered values of the series fi G [0, 1] with precision 5 = 2~ k . All of 
the states of the output register are initially set to zero, so the state of the entire ensemble 
(input and output registers) is given by ® p^ut, 

1 N 

Pin ® Pout = TV E I* >n< An ® |0 >fc< 0|*. (5) 
iV i=l 

Step 2 - Function Evaluation 

The function /, analogous to the oracle in Grover's search algorithm, is evaluated by 
applying a reversible unitary transformation Uf. The transformation has no effect on the 
eigenstates \i > n , but partitions the output register into a set of subensembles, \fo >k , 
defined as the sets of k spins which share the same state. Therefore we have, 

U f \i > n ®|0 > k ^ \i > n ®\fi > k . (6) 
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In a physical implementation, Uj would be the product of a sequence of fundamental unitary- 
transformations for each of the k spins in the output register. The unitary extension of Uf 
to arbitrary output states uses the bitwise XOR operator ©: 

U f \i > n ®|j >fc^ \i >n ®|j © fx >k ■ (7) 

The transformation Uf is applied to the system, in order to evaluate the function / 
simultaneously on the linear combination of all sample points i given by the initial state in 
Eq.(||). This is equivalent to evaluating the function / concurrently on N classical processors 
in a parallel computer. 

Since the initial mixed state is a density operator, the action of Uf can be written as 

1 N 

Uf(pt } ® Pol)U} = - \i >n< An ® I/, >k< Silk (8) 

i=l 

This operation transforms the state of the output register to a mixture that represents 
all of the approximate function values /j, for i = 1, 2, . . . , N. We remind the reader that 
the finite-precision values Si are represented by a binary string of k spins. The states of the 
output register can be written in terms of individual spins, 

\fi >k= \b fi i > ®\b fa > ® . . . ® \b fik >, (9) 

where (6/^, 6/^, • • • , 6/<fe) £ {0,1} are the digits of the approximate function value /j in 
binary format. 

We use the following binary encoding scheme to approximate the set of function values 
S(i) using the 2 k states available in the output register: 



|0i > <g>|0 2 > ® . . 


. ® |o fc 


> 


- /(<) 


e [0,5) 


(10) 


|li > ®|o 2 > ® .. 


• ® |o fc 


> 


- /(<) 


G [5,25) 


(11) 


Id > > ® .. 


. ® |o fc 


> 


- /(<) 


G [25, 35) 


(12) 


1 1 ! > ®|1 2 > ® • • 


• ® life 


> 


- /(<) 


G [1-5,1]. 


(13) 



Alternatively, the approximate function values Si can be defined directly in terms of the 
individual spin values, for example, 

s=5j:y-% v (M) 

3=1 
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to set fi equal to the start of the range intervals given in the encoding scheme above. 

Up to this point, the query complexity of the summing algorithm is one i.e., only one 
function invocation is required. An important, but separate issue that naturally arises 
is whether Uf can be implemented efficiently Despite a query complexity of 0(1), the 
computational complexity of the algorithm can be much higher, if the evaluation of the 
function / is costly. We will not discuss the computational complexity in this paper, since 
we want to keep / as general as possible, but we note that if / is a classically efficiently 
computable function, then Uf can be implemented with comparable complexity, as discussed 



by Nielsen and Chuang []16| . Moreover, functions that cannot be computed efficiently by 
classical devices may be rendered tractable in the future by other quantum algorithms. 



Step 3 - Measurement 

In the last step of the algorithm, we measure the average value of the output register 
in the final ensemble given by Eq.(|8|). This measurement is an analog process, which per- 
forms a single concurrent evaluation, unlike the recursive evaluation on conventional digital 
computers. The result of the measurement is an ensemble average of all the approximate 
function values fi, 

1 N 

/, = ttE/- (15) 

which depends on the distribution of /j and the precision 5 = 2~ k of the encoding. Finally 
the desired sum Sjv,fc can be obtained by multiplying the average value and the total 
number of sample points N. 

Note that in the statement of the problem we have assumed that the number of function 
samples is a power of two, N = 2 n . This ensures that memory resources are employed 
optimally, by using every possible state of the input register to encode the sample points i. 
However, in general, the number of function samples can be arbitrary, in which case only a 
subset of the input register states is used to represent the sample points. 

In the physical implementation of the algorithm, each of the spins in the output register 
generates an output signal jj, proportional to the number of spin subensembles that have 
the j-th spin in the state |1 >. Each signal jj can be transformed into a fraction jj e [0, 1], 



by calibration against the maximum output signal Tj, which is obtained when spin j of the 
output register is set to |1 > for all subensembles. The normalized output signals jj are 
then multiplied by the corresponding binary weight 2 J ~ 1 , cf., Eq.([14]), for j = 1, 2, ... k, to 
give the ensemble average of the output register, 

1 k 

/*= osE 2 *" 1 ^ (17) 

Z 3=1 

and hence the sum Sjv.fc- 

If the measurement sensitivity of the experiment is adequate to distinguish between dis- 
tinct normalized output signals with a precision equal to or better than 1/N, then the query 
complexity remains 0(1). 

However, as the number of sample points N increases, the sensitivity will eventually 
become inadequate. At that point, we cannot reliably measure the output signals for each 
spin in the output register, and significant differences between normalized output signals, 
differences larger than 1/N, will not be detectable in a single experimental trial. We note that 
the error in the ensemble average value fi is given by the weighted sum of the measurement 
errors for each of the spins in the output register, using the exponentially increasing weights 
in Eq.fll7l). To enhance the measurement sensitivity, the algorithm is repeated a number 
of times. For example, in an NMR implementation the proposed algorithm will have to be 
repeated N 2 times, taking into account the square-root scaling of the signal-to-noise ratio 
S oc yN~ e with the number of experimental trials, N e . 



III. DISCUSSION 



The ensemble summing algorithm proposed in this paper uses a radically different ap- 
proach from previous summing algorithms, which rely either on the extrinsic parallelism of 
conventional digital computers, or on the intrinsic parallelism of entangled quantum states. 
Instead, the proposed algorithm uses the parallelism of mixed states in an ensemble of spins 
to evaluate a given function once, and then extracts the measurement result by ensemble 
averaging. We note that the lower bound derived by Nayak and Wu [17|] for the query 
complexity of quantum algorithms that calculate the mean of a function, does not apply to 
our ensemble algorithm, since the polynomial method |18| used in their derivation applies to 
pure states, not mixed states as in the ensemble case. The appropriate lower bound for the 
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(pseudo-pure state) 


(pure state) 


measurement sensitivity scaling 


1/N 


1/N 


1/N 


1 


no. of NMR experimental trials 


N 2 


N 2 


N 2 
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query complexity (single-run) 


0(1) 


O(logiV) 


o(Vn) 


0(VW) 


query complexity (overall) 


0{N 2 ) 


0(N 2 logiV) 


0{N 2 VN) 


0(VW) 



TABLE I: Comparison of measurement sensitivity scaling, number of experimental trials required 
in an NMR implementation, single-run query complexity, and overall query complexity, in the case 
of inadequate measurement sensitivity i.e., for large N. 

ensemble summing algorithm is obtained by considering the query complexity for a classical 
parallel computer with N processors, which is 0(1). 

Table | shows the query complexities of the present ensemble summing algorithm and 
the ensemble search algorithm proposed by Bruschweiler [|R|. We compare both ensemble 
algorithms with Grover's search algorithm || , which provides the critical speedup in existing 
quantum summing algorithms. 

The ensemble summing algorithm has an exponential advantage in terms of query com- 
plexity, relative to the implementation of the quantum summing algorithm using Grover's 
search algorithm with pseudopure states in NMR, as the two algorithms have the same scal- 
ing for the measurement sensitivity. Similarly, the query complexity of the ensemble search 
algorithm is exponentially smaller than that required for Grover's search algorithm using 
pseudopure states in NMR. 

However, a comparison with the (theoretical) implementation of Grover's search algo- 
rithm using pure states shows that both the ensemble summing algorithm and the ensemble 
search algorithm will be more efficient only for a total number of samples below a thresh- 
old value determined by the measurement sensitivity. Bruschweiler |13[ estimates that the 
ensemble search algorithm is more efficient for databases of size N which fulfill the condition 



NVNlog 2 N < S 1 



(18) 



where S is the signal-to-noise ratio of measurements in an NMR implementation. 

The same reasoning leads to the conclusion that the ensemble summing algorithm is more 
efficient than the quantum summing algorithm using Grover's search algorithm with pure 



states, for a number of function samples N given by 

nVn < s 2 , 



(19) 



with respect to the signal-to-noise ratio S. The best available signal-to-noise ratio in present 
NMR technology is S ~ 10 4 , which results in an efficiency threshold value N max k2x 10 5 . 
However for values of N > N max , the query complexity of the ensemble summing algorithm 
is 0(N 2 ), which is greater than the query complexity for quantum and classical summing 
algorithms. 

A different type of error occurs in the initialization step of the algorithm, if a room- 
temperature thermal state is used instead of the initial state p$ in Eq. (|). For an NMR 



implementation, Gershenfeld and Chuang []15| give the thermal state for n spins in the form 

1 N a N 

Pth = T7 X) I* >«< *L + T7 I] Xi\i >n< i\n, (20) 
iv i=l iv i=l 

where a = (« 10~ 6 at room temperature) is the Boltzmann factor of the deviation 

density matrix, and the coefficients \% £ [~ n i n ] represent the net integer sum of various 
spin-up and spin-down combinations of the n spins. 

(n) 

The effect of using a thermal state instead of the equally-weighted mixed state p in be- 
comes present in the measurement step of the algorithm, which results in an ensemble 
average 

1 N a N 

fi = M^fi + (21) 
ly i=l iv i=\ 

The error term in this average imposes an additional limit on the accuracy of the algorithm. 
In the worst case, the difference between the measured and desired values is bounded by 

\fl-U\<™> (22) 

since < n and f\ < 1. In general, the error is much smaller than as the deviation 
matrix is traceless, i.e., J^Xi = 0. The Boltzmann factor a can also be reduced by either 
raising the temperature T, or by choosing a lower average resonant frequency u for the spins. 

The ensemble summing algorithm presented in this paper can be applied to estimating 
the mean and/or the definite integral of a mult i- dimensional function. The validity of 
the algorithm holds for rather general classes of functions, ranging from continuous to the 
Lebesgue measurable and integrable classes, L p , 1 < q < oo. The former is based on the 
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convergence of the Riemann sums to the Riemann integral of "ordinary" functions, while 
the latter is based on the density of simple (Boolean) functions in L p , 1 < q < oo flI5| . 
Unfortunately, these results simply state that in the limit of infinite number of terms, the 
approximating sums coincide with the desired integrals. An evaluation of the error made 
when using a finite number of terms in the sum is impossible, in general. To estimate 
this error, one has to resort to the specific properties of the approximated (sub-classes of) 
functions. For instance, for any Lipschitz function / with Lipschitz constant L, integrated 
over the finite interval [a, b], the error, E^, between the integral of / and the Riemann 
sum with N terms evaluated at equidistant points is bounded by (b — a)L/N. This error 
can now be combined with the error made when estimating the discrete sum with iV terms 
(see Section [II A| ) and an efficient algorithm can then be devised for the estimation of the 
integral. For general functions, the expression of the error as a function of N is unknown, 
although it is known that lim^^ooE^ = 0. However, if decreases with N in a much 
slower fashion, say like En ~ (InN)' 1 , this would translate into a significant increase of 
the number of terms in the sum and therefore an increased complexity, to achieve a given 
overall precision for the integral. The relationship between various functional classes and 
their approximants by Boolean functions is an active research topic that addresses such 
notions as the complexity, capacity, and entropy of a function, which go beyond the scope 
of the present paper. The interested reader is referred to Refs. M, |5], |2D|, |2Tl |22| . 

We conclude by pointing out the two main advantages of the proposed ensemble summing 
algorithm. First, there is no need to maintain quantum coherences for ensemble algorithms, 
so they are easier to implement than their quantum counterparts. Indeed, the scaling of the 
resources required to maintain entanglement in pure state-based algorithms for the duration 
of the computation remains an open question, and could have a potentially large negative 
effect on the exact threshold value N max of the number of function samples. 

Second, for a restricted number of function samples, N <C N max , which is determined by 
the measurement sensitivity, ensemble algorithms may give an exponential speedup over all 
known quantum and classical summing algorithms. In this regime, the proposed ensemble 
summing algorithm requires only a single invocation of the function /. 
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